Antiferromagnetic states and phase separation in doped AA-stacked graphene bilayers 
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We study electronic properties of AA-stacked graphene bilayers. In the single-particle approxima- 
tion such a system has one electron band and one hole band crossing the Fermi level. If the bilayer 
is undoped, the Fermi surfaces of these bands coincide. Such a band structure is unstable with 
respect to a set of spontaneous symmetry violations. Specifically, strong on-site Coulomb repulsion 
stabilizes antiferromagnetic order. At small doping and low temperatures, the homogeneous phase 
is unstable, and experiences phase separation into an undoped antiferromagnetic insulator and a 
metal. The metallic phase can be either antiferromagnetic (commensurate or incommensurate) or 
paramagnetic depending on the system parameters. We derive the phase diagram of the system on 
the doping-temperature plane and find that, under certain conditions, the transition from param- 
agnetic to antiferromagnetic phase may demonstrate re-entrance. When disorder is present, phase 
separation could manifest itself as a percolative insulator-metal transition driven by doping. 
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PACS numbers: 73.22.Pr, 73.22.Gk, 73. 21. Ac 



I. INTRODUCTION 

Graphene is the first experimentally-realizable stable 
true atomic monolayer. It has a host of unusual electronic 
properties 1 - - — . After its discovery, physical properties of 
graphene became the subject of intense scientific efforts. 
In addition to single-layer graphene, bilayer graphene is 
also actively studied. This interest is driven by the desire 
to extend the family of graphene-like materials and to 
create materials with a gap in the electronic spectrum, 
which could be of interest for applications. 

Bilayer graphene exists in two stacking modifications. 
The most common is the so-called Bernal, or AB, stack- 
ing of graphene bilayers (AB-BLG). In such a stacking, 
half of the carbon atoms in the top layer are located 
above the hexagon centers in the lower layer; and half 
of the atoms in the top layer lie above the atoms in the 
lower layer. A different layer arrangement, in which car- 
bon atoms in the upper layer are located on top of the 
equivalent atoms of the bottom layer, is referred to as 
AA-stacked graphene bilayer (AA-BLG), Fig. [TJ So far, 
the most efforts have been focused on studying the AB- 
BLG 4 -, for which high-quality samples are available^. In 
recent years, the experimental realization of the AA-BLG 
has been also reported^—. However, AA-BLG received 
limited amount of theoretical attention^— . 

The tight-binding analysis shows that both AA and 
AB-BLGs have four bands (two hole bands and two elec- 
tron bands). However, the structure of these bands is 
different. In the undoped AB-BLG, two bands (one 
hole band and one electron band) touch each other at 
two Fermi points, and the low-energy band dispersion 
is nearly parabolio 1 ^. The AA-BLG has two bands 
near the Fermi energy, one electron-like and one hole- 
lik o 10 ! 11 . The low-energy dispersion in the AA-BLG is 
linear, similar to the monolayer graphene. Unlike the 



latter, however, AA-BLG have Fermi surfaces instead of 
Fermi points. 



An important feature of the AA-BLG is that the hole 
and electron Fermi surfaces coincide in the undoped ma- 
terial. It was shown in Ref. UH that these degenerate 
Fermi surfaces are unstable when an arbitrarily weak 
electron interaction is present, and the bilayer becomes 
an antiferromagnetic (AFM) insulator with a finite elec- 
tron gap. This electronic instability is strongest when the 
bands cross at the Fermi level. Doping shifts the Fermi 
level and suppresses the AFM instability^. Assuming a 
homogeneous ground state, here we demonstrate that the 
AFM gap A decreases when the doping x grows, and van- 
ishes for dopings above some critical value x c . However, 
the homogeneously-doped state, depending on tempera- 
ture and doping, may become unstable with respect to 
phase separation into an undoped AFM insulator and a 
doped metali£. In the phase-separated state, the con- 
centration of the AFM insulator decreases when doping 
increases. Above a certain threshold value of doping x* , 
the insulator-to-mctal transition occurs. 



In this paper we present a detailed study of the elec- 
tronic properties of the AA-BLG. In Sec. UH we write 
down its tight-binding model Hamiltonian and briefly 
analyze its properties. In Sec. IIII1 we add the on-site 
Coulomb interaction to the Hamiltonian and derive the 
mean-field equations for the commensurate AFM gap 
for finite doping and temperature. The incommensu- 
rate AFM state is analyzed in Sec. IIVI In Sec. [V] we 
demonstrate that the AA-BLG is unstable with respect 
to phase separation within some doping and temperature 
range. The obtained results are discussed in Sec. IVI1 
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FIG. 1: (Color online) Crystal structure of the AA-stacked 
bilayer graphene. The circles denote carbon atoms in the A 
(red) and B (blue) sublattices in the bottom (1) and top (2) 
layers. The unit cell of the AA-BLG consists of four atoms 
Al, A2, Bl, and B2. Hopping integrals t and to correspond 
to the in-plane and inter-plane nearest-neighbor hopping. 



II. TIGHT-BINDING HAMILTONIAN 

The crystal structure of the AA-BLG is shown in 
Fig. [T] The AA-BLG consists of two graphene layers, 
1 and 2. Each carbon atom of the upper layer is located 
above the corresponding atom of the lower layer. Each 
layer consists of two triangular sublattices A and B. The 
elementary unit cell of the AA-BLG contains four carbon 
atoms Al, A2, Bl, and B2. 

We write the single-particle Hamiltonian of the AA- 
BLG in the form 



H = -t { df mA a d mlBa +H, 

(nm)i(T 



(1) 



Here c^ aiaa and d niaa are the creation and annihilation 
operators of an electron with spin projection a in the 
layer i = 1, 2 on the sublattice a = A, B at the position 
n, n is the chemical potential, and (...) denotes nearest- 
neighbor pair. The amplitude t (to) in Eq. ([I]) describes 
the in-plane (inter-plane) nearest-neighbor hopping. For 
calculations, we will use the values of the hopping inte- 
grals t w 2.57 eV, t ~ 0.36 eV computed by DFT for 
multilayer AA systems in Ref. [T3- 

If we perform the unitary transformation 

l _ d nlacr d Tii2aa _ d nlat7 ~ d Tii2aa /-r>\ 

then, Eq. ([T]) can be rewritten as 

H0 = -t^2 (hl Aa h ml3a + H.C.J - (M + t )J2 h laa h naa 
(nm)tr naa 

(nm)(T 



1 (fLlafmBa + H - c ) ~ (A 4 ~ sL3na<7 ■ ( 3 ) 



Therefore, in this representation the Hamiltonian Ho is a 
sum of two single-layered graphene Hamiltoniansi, with 
different effective chemical potential fi ± to ■ 

The Hamiltonian ([3]) can be readily diagonalized. To 
perform the diagonalization, we switch to the fermion 
operators /i kao . and g kac7 , which arc defined in the mo- 
mentum representation, and make the unitary transfor- 
mation 



7kla = 
7k3<r 

where ip^ 



l kAo 



9kAc 



V2 

±9ksA 



V2 
arg (/ k ) 



/k = 1 + 2 exp 



7 k2CT 



7k4<r 



3ik x ao 



h — h p % V* 
rL kA<y ' L kBcr e 



9kAo 



V2 



V2 



cos 



(V3k. 



,(4) 



(5) 



and ao is the in-plane carbon-carbon distance. As a re- 
sult, Hamiltonian becomes 



Ho = 



kstr 



(*) A t 

£ 0k — A* ) iLvlksa 



(6) 



In this equation, the band index s runs from 1 to 4, and 
the band spectra are 



-(i) 

'Ok 

-(3) 
'Ok 



+to 



Kk 



e {2) 

fc 0k 
(4) 



Ok 



+*0 



Kk 

tCk 



(7) 



where £k = |/k|- The band structure obtained is shown 
in Fig. [2] The bands s — 2 and s = 3 cross the Fermi 
level near the Dirac points K = 27r(v / 3, l)/(3v / 3ao) and 
K' = 2tt(V3, -1)/(3\/3oo) [see Fig.^b)]. 

For undoped systems (/.t = 0, half-filling) the Fermi 
surfaces are given by the equation |/k| = to/t. Since 
to ft <C 1, we can expand the function |/k| near the Dirac 
points and find that the Fermi surface consists of two 
circles with radius k r = 2to/(3tao). These Fermi sur- 
faces transform into four circles in doped AA-BLG [see 
Fig. [He)]. 

The most important feature of this tight-binding band 
structure is that at half-filling the Fermi surfaces of both 
bands coincide. That is, the electron and hole compo- 
nents of the Fermi surface are perfectly nested. This 
property of the Fermi surfaces is quite stable against 
changes in the tight-binding Hamiltonian. It survives 
even if longer-range hoppings are taken into account, or 
a system with two non-equivalent layers is considered 
(e.g., similar to the single-side hydrogenated graphene^). 
However, the electron interactions can destabilize such a 
degenerate spectrum, generating a gapi£. 



III. COMMENSURATE 
ANTIFERROMAGNETIC STATE 

The single-electron spectrum described in the previ- 
ous section changes qualitatively when interaction is in- 
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FIG. 2: (Color online) (a) The single-particle band structure 
of the AA-stacked bilayer graphene. It consists of two single- 
layered graphene spectra shifted relative to each other by the 
energy 2to- (b) The k-dependence of the spectra near the 
Dirac point K, located at momentum K. Here, k — K. + Sk y e y . 
The intersection of the bands s — 2 and s — 3 occurs exactly 
at zero energy, which corresponds to the Fermi level of the 
undoped system, (c) The first Brillouin zone (hexagon) and 
the reciprocal lattice unit cell (rhombus) of the AA-BLG. The 
circles around K and K' points correspond to Fermi surfaces 
of the doped system. 



order parameters as 



A 



1A 



A, 

A2B = — Aie 



±2A 



A. 



(9) 
(10) 



and A is real. Such AFM order, when spin at any given 
site is antiparallel to spins at all four nearest-neighbor 
sites, is called in the literature G-type AFM. Other types 
of spin order are either unstable or metastablc. 

In the mean-field approximation, the interaction 
Hamiltonian has the form 



TJ- 



MF 



jV 



4A 2 



U(n 2 



1 



Ux x -v 



where x = n — 1 is the doping level, n is the number of 
electrons per site, and Af is the number of unit cells in 
the sample (a unit cell of AA-BLG consists of four car- 
bon atoms, see Fig. [T]). Below, when quoting numerical 
estimates for doping, we will write percentage of 

the total number of carbon atoms in the sample. 
We introduce the four-component spinor 



lb 



Aa ' AtlAcs ' ^klB<r i d \a1Ba ) ' 



(12) 



which can be used to build an eight-component spinor 
= Met) V'kJ.) - ^ n t crms °f this spinor, the mean field 



Hamiltonian H MF = H 



H^F can be written as 



eluded. Specifically, using mean field theory, we will 
demonstrate that the degenerate Fermi surface of the un- 
doped AA-BLG is unstable with respect to the sponta- 
neous generation of AFM order. 



A. Mean-field equations 

We approximate the electron-electron interaction by 
the Hubbard-likc interaction Hamiltonian: 



where 

4A 2 JJr 
E = ^-U(n*-1), // = /;-^. (14) 

In these equations, Eq is a c-numbcr, fjf is the renormal- 
ized chemical potential, -f/ok, and A are 4x4 matrices 



Hi 



Ok 



t tf u \ 

to tf k 

tf* t 

*/; t / 



(15) 



H 



int 



4e 



(8) 



where n niaa = dnj a(7 d nia(T , and a = -a. It is known 
that the on-site Coulomb interaction U in graphene and 
other carbon systems is rather strong, but the estimates 
available in the literature vary considerabl y 20 ' 21 , ranging 
from 4-5 to 9-10 eV. 

We analyze the properties of the Hamiltonian H = 
Ho + -ffint in the mean-field approximation. We choose 
the x-axis as the spin quantization axis, and write the 



A = 



/-A 

A 

A 

V 



(16) 



We diagonalize the 8x8 matrix in Eq. (fT3")) and obtain 
four doubly-degenerate bands 



,(1,4) 



= TJA2 + (t&+toY 



.(2,3) 



WA 2 + (tCk-tor 



(17) 
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To determine the AFM gap A we should minimize the 
grand potential Q. The grand potential per unit cell is 



rfk 



In 



,(m'-4 s) )/ t 



(18) 



where Vbz is the volume of the first Brillouin zone. 

To evaluate integrals over the Brillouin zone it is con- 
venient to introduce the density of states 



(19) 



This function is non-zero only for < C < 3. It 
is related to the graphene density of states p gr {E) as 
Pgr (E)=p (\E/t\)/t (seeRef.H). 

Minimization of D, with respect to A gives the equation 



U 
4t 



d( P o((){F(VS 2 + (( + (o) 2 ) + 



F (V5 2 + (C - Co) 2 ) 



(20) 



where S = A/f, Co = to/t, and 



F(e) 



f(-te - p') - f(ts - p!) 



, f(E) 



E 

eT 



1 



(21) 

Equation (|20|) determines the gap A as a function of the 
renormalizcd chemical potential p' . To find A as a func- 
tion of doping, we need to relate the doping and the 
chemical potential. It is easy to prove that 



1 d(Q - E ) 
n=1 + X = -4^p^- 

Then, using Eqs. (fT5j) and (fTO|) we derive 



x = \ Jd(po(C)[G (V<S 2 + (C + Co) 2 ) + 



(22) 



G (V$ 2 + (C- Co) 2 
where G(e) = f(-ts - //) + f(te - p!) - 1 



(23) 
(24) 



Solving Eqs. ([20]) and (J23J) we obtain the AFM gap 
A(x,T) and the chemical potential p(x,T). 

The solutions of Eqs. ((20]) and ((23]) satisfy the fol- 
lowing relations: A(— x,T) ~ A(x,T) and p(—x,T) — 
—p(x,T). They are consequences of the particle-hole 
symmetry of the model Hamiltonian. The next-nearest- 
neighbor hopping breaks this symmetry. However, our 
analysis shows that corrections, introduced by these 
terms, do not exceed 1-2% for the range of parameters 
characteristic of graphene systems. Assuming particle- 
hole symmetry, below we only consider electron doping, 
x > 0. 




£/[eV] 



FIG. 3: (Color online) (a) The AFM gap ratio A/A ver- 
sus doping x/x c for different values of the on-site Coulomb 
repulsion U : (red) squares correspond to U — 5.5 eV, (blue) 
circles to U = 7eV, (green) triangles to U = 9eV. The solid 
(black) curve is A(x)/A = y/l - x/x c . (b) The dependen- 
cies A(x = 0, T = 0) = Ao and critical doping x c versus U. 
The solid curves are numerical solutions of Eqs. (|25|l and (|26|) . 



while the dashed curves are calculated using the approximate 
analytical solution, Eqs. (|29p and (|30p . 



B. Zero temperature 
If T = 0, Eqs. (20) and (23) become 



i - e L'/t - yW(c + Co) 2 ) 



47 / 



l-e [p'/t- J S 2 + ((-{ ) 



+ 



6 2 + ((+( Q y 



VW(C-Co) 2 



(25) 



© (Va - \A 2 + (C + Co) 2 ) + 



e(yA-\A 2 + (C-Co) 2 ) 



(26) 
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where Q(x) is the step function. At half-filling, n = 1, 
x = 0, and /i = // =0. The lower two bands are filled, 
the upper two bands are empty, and both 6-functions in 
Eq. (2"5)l are zero for any (. 

When doping is introduced, analysis of the latter 
equations shows that /z' changes abruptly from zero to 
the value \rf > A. The gap A(x,T = 0) decreases 
monotonously from A(x = 0, T = 0) = A to 0, when 
x increases from to some critical doping x c . To find 
x c we must put 5 = into Eqs. (|25|) and (|26j) and solve 
them for fi' and x = x c . 

Equations (f25|) and (f26| can be solved analytically, if 
Ao <C t,to. Using the asymptotic expansions of integrals 
in these equations for small 5 we obtain^ 



A(x,0) 
/i(x,0) 



Ao, 1 



■ ' r. 



A 
A po(Co 



Sgn(x) -2^ 



A i 



(27) 
(28) 

when t < t . (29) 



2t 7rV3i 2 
In this limit the value of A is given by the relation-^ 

4t~U V (( ' 



A = 2^/t (3t - t )exp 
where 



2(7po(Co) 



v(Co)= d{ 



Po(C) . Po(C) - A) (Co) 



C + Co 



IC-Co 



(30) 



(31) 



The dependence of the ratio A(x,0)/Ao on x/x c for 
different values of U is shown in Fig. (3{a) . Figure [3Jb) 
shows A and functions of U calculated both nu- 

merically [Eqs. (55) and (H] and analytically [Eqs. (M]) 
and pO])]. Equations (J28J) , j29]) together with Eq. <J30]) 
for Ao are valid if C/ < 6 eV. However, Eq. ((2T)) is accu- 
rate for any U, provided that x c and A are calculated 
numerically from Eqs. (|23t and (|2^|) [sec Fig. EJa)]. 



C. Finite temperatures 

In this subsection we will analyze the finite- 
temperature solutions of the mean field equations (f2"0)) 
and (|23[) . However, it is necessary to remember that in 
2D systems no long-range order is possible if T > 0. In 
such a situation the mean field solutions characterize the 
short-range order, which survives for sufficiently low T. 
The effects beyond the mean field approximation will be 
discussed in subsection IIII D[ 

Solving numerically the mean field equations (|2H)) 
and (|23p , we find A as a function of doping x and temper- 
ature T (see Fig. H]). The temperature Tmf at which A 
vanishes is the mean field transition temperature (see in- 
set of Fig. U). The transition temperature, as a function 




FIG. 4: (Color online) The dependence of A(x,T) on doping 
x calculated for U = 5.5 eV and different T/A : (1) T/A = 
0.06, (2) T/Ao = 0.17, (3) T/A = 0.33, (4) T/A = 0.41, 
(5) T/A = 0.47, (6) T/A = 0.52, (7) T/A = 0.55, and 
(8) T/Ao = 0.58. Inset: The dependence of the mean-field 
transition temperature on doping. The re-entrance from PM 
to AFM state exists in the doping range x c < x < 1.231:Cc; 
x c = 0.128% and A = 0.124 eV. 



of doping x, is not a single-valued function. Instead, it 
demonstrates a pronounced re-entrant behavior. We dis- 
cuss this unusual phenomenon in more details in Sections 
13 El andEH 

Equations (f2T)]) and (|2"3"]) can be simplified in the case of 
small gap, when Ao <C to, t. Neglecting terms of the or- 
der of Ag/i 2 in Eq. (J2DJ and taking into account Eq. ([50)) 
for Ao, we obtain the following equation for A 



In 



A 



1 



dz acosh 



zT 
A" 



(32) 



A/T 



cosh 



H'/T 



cosh 



z + n'/T 



In the same limit, we derive from Eq. (|23|) the relation 
between jj,' and x in the form 



x 



T 



, dz\l Z 2 

2A 0l / V T 2 



(33) 



cosh 



A/T 

-2 l'z ~ fJt' '/T 



cosh 



z + fi'/T 



At half-filling (x = 0) we find from Eq. (32]) the BCS- 
like result Tmf(0) = 0.567Ao. If we normalize x by x c 
and A, fj/ and T by A , then, Eqs. (22) and O do 
not include any parameter characterizing the AA-BLG 
band structure. Thus, if the electron interaction U is 
not large, the obtained results do not depend on details 
specific to the AA-BLG and are valid for other systems 
with imperfect nestin g 19 ' 22 ' 23 . 
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D. Crossover temperature 

In 2D systems, finite-temperature fluctuations destroy 
the AFM long-range order. Then, the results obtained 
above in the mean-field approximation are valid only if 
the mean-field correlation length £ = Vf/A is smaller 
than the spin-wave correlation length £ su , (here vp = 
3aoi/2 is the Fermi velocity in our model). Otherwise, 
short-range ordering of spins disappears, and we cannot 
define the AFM order even locally. 

In the limit £ sw > £, the spin fluctuations can be de- 
scribed using the nonlinear tr-model with Lagrangia n 2 ' 25 



£sw = §[(d t D) 2 -c 2 w (d r D) 2 



(34) 



where D is the unit vector along the local AFM magne- 
tization. The spin-wave stiffness p and velocity c sw can 
be evaluated from Eqs. (7.89) and (7.90) of Ref. Eg 



V2' 



= fto/(87raf), 
P |A/(16n 



if t > A, 
if t < A. 



(35) 



The correlation function K(r) = (D(r)D(O)) can be ob- 
tained using the Lagrangian (j34|) . At large distances it 
behaves a s 24 ' 25 



K(r) « 1 



T 



In 



3aot 



(36) 



where 7 is the Euler's constant. The spin-wave correla- 
tion length £ sw describing the characteristic size of the 
short-range AFM order can be estimated using the equa- 
tion A'(£ sw ) = 0. Thus, we have 



6 



a t 



(37) 



Solving the equation £ sw = £, we find the crossover tem- 
perature T* between the short-range AFM and the PM. 
The short-range AFM order exists over distances of about 
Csw > do, if T < T*, and it is destroyed if T > T*. Our 
numerical analysis shows that T*(x)/Tmf(x) ~ 0.8-0.9 
for any ratio A/t. Thus, the mean-field transition tem- 
perature gives an appropriate estimate for the AFM-to- 
PM crossover temperature. 



IV. INCOMMENSURATE 
ANTIFERROMAGNETIC STATE 

The G-type AFM state considered above has the small- 
est value of the grand thermodynamic potential Q among 
other states with commensurate magnetic order. How- 
ever, further optimization of fl could be achieved if we al- 
low the local direction of the AFM magnetization slightly 
rotate from site to site 1 ^. Then, the translation invari- 
ance with a lattice period disappears. Such a state is 




x/x 1.0 



FIG. 5: (Color online) The dependence of A (red solid 
curve) and |q| (blue dashed curve) on doping x calculated 
for T/Ao = 0.06. The dot-dashed curve is the gap A calcu- 
lated for q = and U — 8eV. The doping x is normalized 
by the critical doping x c , calculated for the commensurate 
AFM state. The incommensurate AFM exists in a slightly 
larger doping range than the commensurate AFM. Note that 
the gap for commensurate AFM remains non-zero even for 
x > x c . This is a manifestation of the re-entrance, see inset 
of Fig. H 



referred to as incommensurate (or helical) AFM. The 
complex order parameter for this state has a form 



U(d 



i f d 



p* qn A 



(38) 



where q describes the spatial variation of the AFM 
magnetization direction, the position vector n specifics 
the location of a given carbon atom, and Ai a satisfies 
Eqs. (fl"0|) . The averaged electron spin S n?a at site n lies 
in the xy-plane. It is related to the order parameter as 

( d Lat d mai) = 5 ma + ^nia' As a r6Sult > We ° btain 



U 



(cos(qn), sin(qn)) 



(39) 



The mean-field version of the interaction Hamilto- 
nian ((5J corresponding to the order parameter A n i a , 
Eq. (|38j) . can be written in the momentum representa- 
tion as [c.f. Eq. (fTT])] 



4A 2 



V 



C/(n 2 



Ux 



kz« 



E 

kiaer 



Tiki, 



(40) 



k+^iat" k— Sia4- 



:— -^ia\. k+^?'af 

It is convenient to redefine the spinor ^ (sec Sec. IIII Al) : 

iHi^fi^)' ( 41 ) 

We can rewrite the mean- field Hamiltonian Hq + H^jf in 
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the form 



H, 



Ok 



kql 



q 
A 



A 



H, 



Ok-Sf 



kq' 



(42) 

where /f k an d A are given by Eqs. (fT5|) and (|T5|) . re- 
spectively. 

The electron spectrum in the incommensurate AFM 
state is found by diagonalization of the 8x8 matrix in 
Eq. (|4"2"j) . It consists of 8 non-degenerate bands, , 

(s) 

s = 1, 2, . . . , 8. The analytical expression for E^. ' can be 
obtained in the limit |q| <C 1/dOi 



E 



t 



A2 



+ 



^ 0± C k+ ^+C k .q 
t 2 



(43) 

If q = 0, this spectrum coincides with the spectrum of 
Eqs. (H7D- 

The expression for the grand potential Q has similar 
structure as Eq. ()18[) . but now the summation includes 
eight bands: 



1 



(44) 



Minimization of with respect to A and q, together with 
the condition relating x and fx' gives the closed system of 
equations for calculating A(x,T), q(x,T), and /x(x,T): 



5A 



0. 



9f2 
9q 



0, 1 



a(n - g ) 

dfi' 



(45) 



We calculate the functions A(x,T), q(x, T), and 
/Lt(x, T) numerically for different values of U. Typical 
curves A(x) and |q(x)| are shown in Fig. [S] For com- 
parison, the curve A(x) calculated for the commensu- 
rate AFM is also plotted. We see that the incommensu- 
rate AFM state exists in a slightly wider doping range 
than the commensurate one. The incommensurate phase 
arises at arbitrary small doping if T = 0. At non-zero 
T the commensurate AFM state is stable until dop- 
ing exceeds some T-dependent threshold x q (T). The 
curve T q (x) separates the incommensurate and commen- 
surate AFM states; the more symmetrical AFM state 
with q = lies above T q (x). 

The phase diagrams of the model in the x-T plane 
are shown in Fig. [6] for three different values of U. The 
diagrams for small U (< 6eV) and large U (> 6eV) 
demonstrate a qualitative difference. Namely, for small 
U [Fig.[B{a), U = 5.5 eV] the re-entrance, seen in the in- 
set of Fig. |4l disappears. It is masked by the incommen- 
surate AFM phase. For larger U, however, it survives, 
Fig. H)Ib,c). Re-entrance is an unusual phenomenon be- 
cause the ordering occurs as the temperature increases. 
If re-entrance is a genuine feature of the model, or it is 
an artifact of the mean field approximation, whose relia- 
bility deteriorates when U grows, we do not know. Simi- 
lar behavior was predicted theoretically for quarter-filled 



0.08 - 

> 

u 

E-H 

0.04 



fa) 


1 




1 




L 

l/=5.5eV |- 


AFM 








q=0 




>V\. 


PM 




PS 


\ X. 
\ X 

'AFT 












X 



X[%] 0.2 



0.2 h 
> 

E-H 

0.1 



(b) 


i 1 i , 


1 i 




f/=6.5 eV | 


AFM 






q=0 








AFM // 




/ 

/ 


q * / / 


PM 


/ 

/ / 


\ : / 
i \ , ■ 1 i 


i 



x [%] 2 



0.2 
> 

0.1 



(c) 

'AFM / 
q=0 /' 



U=S eV 



AFM 
q*0 



;ps 



PM 



5 x [o/ ] 10 

FIG. 6: (Color online) The phase diagram of the model in 
the (x, T) plane, calculated for the electron doping x > and 
U = 5.5 eV (a), U = 6.5 eV (b), and U = 8eV (c). Solid 
(red) curves are Tmf(x), (blue) dashed curves are T q (x), at 
which the commensurate-incommensurate transition occurs. 
The dotted (red) curves are Tmf(x), calculated without tak- 
ing into account the incommensurate AFM state. The dot- 
dashed (green) curves show the region of phase separation. 
For hole doping (x < 0) the results are the same. 



Hubbard model at moderate interaction strength^, and 
numerically for classical rotor model^. 
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V. PHASE SEPARATION 

In our discussion above we implicitly assumed that the 
ground state of the AA-BLG is spatially homogeneous. 
However, this is not always true: it was predicted in 
Ref. [l6| that there is a finite doping range where the AA- 
BLG separates in two phases with unequal electron den- 
sities 711,2 = 1 + xi.2- Indeed, if Ao <C t,to we can use 
Eq. ([28"]) and obtain that 



dx 



<0, 



if — < . 

t t 



(46) 



The negative value of the derivative d/j,/dx indicates the 
instability of the homogeneous state toward phase sepa- 
ration^. 

If the possibility of the incommensurate AFM is ig- 
nored, then a zero-temperature phase separation^ oc- 
curs between the AFM insulator (x± = 0) and the PM 
(U < 6cV) or the AFM (U > 6cV) metal (x 2 > 0). 
Here we study phase separation taking into account the 
incommensurate AFM phase and non-zero temperature. 
We numerically analyze the stability of the homogeneous 
state using the dependence of the chemical potential /i 
on the doping x. 

A typical dependence n(x) for non-zero temperature 
is shown in Fig. [7j The derivative dfi/dx is negative in 
some range of doping and the system separates in com- 
mensurate (q = 0, X\ < x) and incommensurate (q ^ 0, 
X2 > x) AFM phases. The doping concentrations x\ 
and X2 are found using the Maxwell construction^: the 
(black) horizontal line is drawn in such a manner that the 
areas of the shaded regions in Fig. [7] are equal to each 
other. When temperature increases, the doping range 
x\ < x < X2, where the phase separation exists, be- 
comes narrower and disappears at some critical temper- 
ature. Our calculations show that the separated phases 
are AFM with q = and q ^ for any values of the 
model parameters. The region of the phase separation in 
the (x, T)-phase diagram is shown in Fig. [5] by (green) 
dot-dashed lines. 



VI. DISCUSSION 

In this paper we study the evolution of the electron 
properties of AA-BLG with doping x and temperature 
T. We calculate the phase diagram of the system in the 
(x, T)-plane. This diagram includes regions of the AFM 
commensurate and incommensurate states, a region of 
phase separation, and the PM state. With good accuracy, 
the electronic properties of the AA-BLG arc symmetric 
with respect to the electron (x > 0) or hole {x < 0) dop- 
ing. The maximum crossover temperature between the 
short-range AFM and the PM states depends on the on- 
site Coulomb repulsion U. For example, AFM ordering 
can exist up to temperatures of about 50 K if U = 5.5 eV 
and to temperatures much higher than room tempera- 
ture if U > 6.5 eV. At present, there is no consensus on 




FIG. 7: (Color online) Chemical potential /j, of the homoge- 
neous state versus doping x; U — 5.5 eV and T = 0.014 eV. 
The vertical dot-dashed line separates the AFM states with 
q = and q ^ 0. In the doping range x\ < x < X2, phase 
separation occurs. The values 2:1,2 are determined by the 
Maxwell construction: the horizontal (black) line is drawn in 
such a manner that the areas of the shaded regions are equal 
to each other. 



the value of the on-site Coulomb repulsion in graphenc- 
bascd materials. However, it is commonly accepted that 
U lies is in the range 6 < U < 10 eV and, consequently, 
the AFM state can be observed in the AA-BLG. 

The critical doping value x c , at which the AFM is re- 
placed by the PM, also strongly depends on U, chang- 
ing from - 0.1% if U = 5.5 eV to ~ 10% if U = 
8eV. For graphene systems, doping of about 10% and 
even higher was achieved 3 -^. Similar to single-layered 
graphene, the AA-BLG can be doped by using appro- 
priate dopant a 30 ' 31 , choosing the substrate and applying 
a gate voltag e - 32 : 33 , or by combinations of these methods. 

We predict the existence of phase separation in the 
AA-BLG. The separated phases have different electron 
concentrations, x\ and X2, and the phase separation will 
be frustrated by long-range Coulomb repulsion 3 ^. In 
this case the formation of nano-scale inhomogeneities is 
more probable. The electron-rich phase (incommensu- 
rate AFM) is metal and the electron-poor phase (com- 
mensurate AFM insulator if T = 0) is insulator or "bad" 
metal. Thus, the percolative insulator-metal transition 
will occur when the doping x exceeds some threshold 
valued, which is about 0. 5(2:1+0:2) in 2D systems. Phase 
separation exists in the doping range x\ < x < X2, and 
^2^1% for any value of U . Depending on U, phase sep- 
aration could be observed from 30-40 K to room and even 
higher temperatures (see Fig. [6]). This makes AA-BLG 
promising for applications. 

The incommensurate AFM phase is mathematically 
equivalent to the Fulde-Ferrel-Larkin-Ovchinnikov state 
in superconductors 3 ^ 3 ^, which is sensitive to disorder— 
and difficult to observe experimentally. Consequently, it 
is reasonable to expect that the incommensurate AFM 
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phase can be destroyed by factors our study did not ac- 
count for. Our calculations predict that in this case the 
region of phase separation changes only slightly in the 
phase diagram. However, the separated phases would be 
the AFM insulator and PM (U < 6eV) of AFM metal 
(U > 6eV). 

To conclude, we studied the phase diagram of the AA- 
stacked graphene bilayers on the doping-tcmpcrature 
plane. It consists of paramagnetic and antiferromagnetic 
(both commensurate and incommensurate) homogeneous 
phases. In addition, a region of phase separation is also 
identified. Magnetic properties of the AA-BLG may sur- 
vive even at room temperature. 
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